TUTORIAL

Mapping Voxel-based niftis to ROI-based surface data

prepare environment and load libraries

Make sure you have installed the latest version of the connectome workbench, which you can download here: https://www.humanconnectome.org/software/connectome-workbench

#### load required libaries ####
require(neurobase)
require(ciftiTools)
require(plyr)
require(dplyr)

#### set path to workbench (example for a MacOS environment) ####
ciftiTools.setOption('wb_path', '/Applications/workbench/bin_macosx64')

## Using this Workbench path: '/Applications/workbench/bin_macosx64/wb_command'.

extract ROI-based data

For this example, we will extract ROI based data from a tau-PET image using the 17 network version of the Schaefer 200 ROI atlas. The atlas including the reference (Schaefer et al., Cerebral Cortex, 2018) can be found here: https://github.com/ThomasYeoLab/CBIG/tree/master/stable_projects/brain_parcellation/Schaefer2018_LocalGlobal

set tutorial path

ste the path to your downloaded tutorial

sdir = "/Users/isd_neuroimaging/Library/Mobile Documents/com~apple~CloudDocs/Github/Tutorials/Connectome_workbench/Map_nifti_to_ROIs/"

read in images

read in atlas

# load the atlas
atlas_readin <- neurobase::readnii(
  file.path(sdir, "/Example_data/rSchaefer2018_200Parcels_17Networks_order.nii.gz")
)
neurobase::ortho2(atlas_readin)

read in the PET image

# load the atlas
pet_readin <- neurobase::readnii(
  file.path(sdir, "/Example_data/smoothed_8mmFWHMwsub-0555_ses-2020-02-06_pet-AV1451_flirtT1reg_INFCER_SUVR.nii.gz")
)

## Warning in datatyper(nim, warn = warn): Need to change bitpix and datatype to
## FLOAT64 due to NAs

neurobase::ortho2(pet_readin)

check overlap between atlas and pet image

# check overlay
neurobase::ortho2(pet_readin, atlas_readin)

# check if dimensions match
dim(pet_readin) == dim(atlas_readin)

## [1] TRUE TRUE TRUE

extract ROI-based mean values

pet_atlas_based = 
  data.frame(pet = as.vector(pet_readin@.Data),
             atlas_idx = as.vector(atlas_readin@.Data)) %>% 
  filter(atlas_idx !=0) %>% 
  group_by(atlas_idx) %>% 
  summarize(mean = mean(pet))
head(pet_atlas_based)

## # A tibble: 6 × 2
##   atlas_idx  mean
##       <int> <dbl>
## 1         1 1.11 
## 2         2 0.742
## 3         3 0.860
## 4         4 0.904
## 5         5 1.09 
## 6         6 1.16

# show distribution of SUVR data
hist(pet_atlas_based$mean)

map to xifti file

# read in atlas
atlas_dscalar <- read_xifti(file.path(sdir, "/Example_data/Schaefer2018_200Parcels_17Networks_order.dscalar.nii"))
atlas_dscalar_vec <- c(as.matrix(atlas_dscalar))

# determine values to map to surface
vector_to_map = pet_atlas_based$mean

# map ROI-based PET data to cifti
xii_to_map <- c(NA, vector_to_map)[atlas_dscalar_vec + 1]
xii1 <- select_xifti(atlas_dscalar, 1)
xii_to_map <- newdata_xifti(xii1, xii_to_map)

# view cifti directly in R
view_cifti_surface(xii_to_map, colors = "magma", zlim = range(vector_to_map),
                   borders = T)

write cifti to disk

ciftiTools::write_cifti(xii_to_map, file.path(sdir, "/Example_surface"))

## Writing left cortex.
## Writing right cortex.
## Creating CIFTI file from separated components.